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Abstract. In this paper we present a new local La- 
grangian approximation to the gravitational dynamics of 
cold matter. We describe the dynamics of a Lagrangian 
fluid element through only one quantity, the deformation 
tensor. We show that this tensor is clearly suited to the 
study of gravitational dynamics and, moreover, that know- 
ing its evolution is enough to completely describe a fluid 
element. In order to determine this evolution, we make 
some physical approximations on the exact dynamical sys- 
tem and we thus obtain a closed local system of differen- 
tial equations governing the evolution of the deformation 
tensor. Our approximate dynamics treat exactly the con- 
servation of mass, of the velocity divergence and of the 
shear and is exact in the case of planar, cylindrical and 
spherical collapses. It also reproduces very accurately the 
evolution of all dynamical quantities for a very wide class 
of initial conditions as illustrated by a detailed compar- 
ison with the homogeneous ellipsoid model. Beside, we 
highlight, for the first time, the important dynamical role 
played by the Newtonian counterpart of the magnetic part 
of the Weyl tensor. 
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1. Introduction 

Two main analytical approaches have been used to 
study non-linear gravitational instability dynamics in 
cosmology. The first one is Eulerian, the funda- 
mental quantities are then the density and veloc- 
ity fields evaluated at the comoving Eulerian co- 
ordinate X (Goroffetal. 1986, Grinstein & Wise 1987). 
The second method is that of Lagrangian trajecto- 
ries (Zeldovich 1970, Moutarde et al. 1991, Buchert 1993, 
Bouchet et al. 1995). However, since the equations gov- 
erning the motion of a self-gravitating dust are highly non- 
linear and non-local, finding an exact solution in the gen- 



eral case, both in the Eulerian or Lagrangian approaches 
is impossible. Consequently, one is then forced either to 
restrict the problem to particular geometries (spherical or 
planar) or to use approximation techniques such as a per- 
turbative theory. In the first case, one can study the very 
nonlinear evolution of a single object, but the range of 
initial conditions is then by definition very limited. The 
second case applies to any type of initial conditions, it con- 
sists of expanding the solutions of dynamical equations. In 
the Eulerian formalism one develops the density contrast 
and peculiar velocity fields as follow: 



n=0 



E 

n=0 



(1) 



These developments are introduced in the dynamical 
equations (Mass conservation, Euler and Poisson equa- 
tions). The equations are then solved keeping only terms 
in e and after in and so on. The Eulerian expansion of 
the density contrast is valid when S < 1, but when 5 « 1, 
neglecting the term of high order is no longer possible. 
In the Lagrangian formalism the solution is obtained as a 
perturbated trajectory about the linear displacement 



(2) 



n=0 



where q is the Lagrangian coordinate. This kind of ex- 
pansion allows one to follow the dynamics longer than in 
the Eulerian description, since the Lagrangian picture is 
intrinsically non-linear in the density field. The first order 
recovers the Zeldovich approximation (Zeldovich 1970). 
Higher orders provide corrections but the calcula- 
tion rapidly becomes very tedious (Moutarde et al. 1991, 
Buchert 1992) and their domain of application are in fact 
restricted to the quasi- linear regime {d less than a few). 

In this paper we propose a new Lagrangian approxi- 
mation. It allows to compute the local properties of a fluid 
element during its motion, eventhough it does not deter- 
mine its trajectory. From this point of view, our work is 
in the same spirit as that of Bertschinger & Jain (1994) , 
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however, in our case, the nonhnear dynamic will be stud- 
ied through only one quantity, the deformation tensor 



dxi 



(3) 



X and q are respectively the Eulerian and Lagrangian co- 
ordinates of a fluid element. 

The deformation tensor when it is not singular, 
contains all the dynamical information on a fluid el- 
ement excepted its trajectory. From this tensor, and 
from its time-derivatives, one can compute the den- 
sity field, the velocity divergence, the shear tensor 
and the tidal tensor. Consequently, the influence of all 
these quantities which are very operative during the 
evolution of a fluid element (Bertschinger & Jain 1994, 
R. van dc Weygacrt & A. Babul 1994), will then be natu- 
rally taken into account in our formalism. Moreover as Wij 
are by definition always smaller than 1 (in absolute value) , 
for a collapsing fluid element, they are very suitable quan- 
tities for all approximate perturbative description and can 
allow one to follow the evolution of large density contrasts 
or large displacements, this is not possible with the usual 
Eulerian or Lagrangian perturbative approach. If a fluid 
element is expanding in one direction the corresponding 
Wij can be greater than 1. In that respect, the deforma- 
tion tensor is more adapted to describe collapsing regions 
than expanding voids. 

We deduce the nonlinear dynamics of the deformation 
tensor by using an alternative procedure to the usual per- 
turbative expansion of the solutions of the exact dynam- 
ical system. Our procedure consists, after some physical 
approximations, in deducing directly from the exact sys- 
tem a new approximate dynamical system where the non- 
locality of the gravitational dynamic introduced by the 
Poisson equation has been cancelled. The new closed set of 
local Lagrangian equations is then solved exactly. It takes 
into account nonlinear effects which could be absent in a 
perturbative solution. In other fields of research similar 
procedures have been used. It is the case for example, in 
the turbulent hydrodynamic or in the state solid physic, 
where some important physical processes, as the gener- 
ation of solitons and the tunnel effect respectively, can 
not be explained by a perturbative solution. Our system 
reproduces exactly the spherical, cylindrical and planar 
collapses, but it also allows one to describe a much wider 
class of initial conditions and to follow the dynamics into 
the highly nonlinear regime. 

A comparative analysis of our approximation with 
the results of Bertschinger & Jain 1994 and Lagrangian 
perturbative solutions has also been performed. It has 
allowed a better understanding of the respective influ- 
ence of different dynamical quantities during the col- 
lapse, as for example the Newtonian limit of mag- 
netic part of Weyl tensor which has been extensively 
discussed in recent papers (Bertschinger & Jain 1994, 



Bertschinger & Hamilton 1994, 

Kofman & Pogosyan 1995, Ellis & Dimsby 1995) 

The outline of the paper is the following: In section 
2 we present the dynamical equations of the deforma- 
tion tensor and our approximations on these equations. 
In section 3 we give a physical explanation to our approx- 
imations and evaluate their accuracy. We also suggest a 
classification of collapses according to their geometrical 
properties. We compare our results with other approxi- 
mate descriptions in section 4 Finally, conclusions and a 
summary of our results are given in section 5 The com- 
plete perturbative solutions to our dynamics is given in 
the appendix. 

2. Dynamics of the Deformation Tensor 

The evolution of the cold matter in an expanding Universe 
is governed by the following set of equations: 



^ + -Vx.(l + <5)v = 
ot a 

1 . „ ^ a 1 „ 

^ + - v.Vx)v +-v = — Vx$ 

ot a a a 



(4) 

(5) 
(6) 



<& is the gravitational potential of the fluctuating part of 
the matter density, (p is the average matter density, S is 
the density contrast) v is the peculiar velocity (v = o^) 
and a is the scale factor. The non-locality introduced by 
the Poisson equation (6) makes this system very difficult 
to integrate. However, by using the deformation tensor 
J, it becomes possible to derive an approximate local de- 
scription. 

First, we change from the Eulerian to the Lagrangian 

time derivative and we introduce the conformal time t 
defined by dr = the previous set of equations can 
then be written as: 

■a{l + S)e = (7) 



dr 



dv 



^ = -Vx$ = +dv + a2(v.V)v 



(8) 
(9) 

where <!> = 0^$ and 6 is the velocity divergence. 

The conservation of mass (equation (7)) can be rewrit- 
ten in terms of JT" as 1 -|- ^ = j where J = det J. The 
derivative with respect to q of equation (8) and the prop- 
erties of commutation between the operators Vq and ^ 
yields 



(10) 



where * denotes a matrix product. Finally, by introduc- 
ing Tij = Vj;.Vj;^$ in equation (9), the new system of 
equations becomes: 

{dlJ = -J*J^ (11) 
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We have defined /3(t) = 47rGa^(r)p(T). All the effects of 
n on the dynamics arc contained in this function /?. In an 
= 1 universe (3 = These equations were first writ- 
ten by Lachize-Rey (Lachieze-Rey 1993) who proposes a 
solution where the deformation tensor has the particular 
form: J{c[,t) = £'q(t)j7o(q)- ^7 is a tensor field which 
depends only on the Lagrangian coordinates and -Dq(t) 
is a scalar function of time which is given by a differen- 
tial equation whose coefficients are functions of the eigen- 
values of j7o(q) • Assuming that the deformation tensor 
remains proportional to its initial value during the evolu- 
tion of the fluid element greatly simplifies the dynamics. 
However, this is a very restrictive assumption, and such a 
solution admits a perturbative expansion which generally 
differs even at second order from the usual Lagrangian 
perturbative solution. 

In the general case has nine independent compo- 
nents and system (11) cannot be solved. However, if the 
flow is irrotational then Sij = is symmetric and 

since S = -J^^^ * d-^J , this implies that J7 is also sym- 
metric, with only six independent components. In the rest 
of this paper we restrict ourself to irrotational flow. More- 
over, if J7 is symmetric, it can be diagonalizcd. Wc per- 
form this diagonalization in the initial conditions and we 
consider that the first and second order time-derivative of 
non-diagonal components are initially vanishing, in order 
that J remains diagonal during all the evolution. This 
assumption restricts somewhat the range of initial condi- 
tions. However, as it has been shown by Moutarde et al. 
(1991) in the case of the evolution of three caustics the 
non diagonal terms which appear during the evolution, do 
not modify strongly the dynamics of the diagonal terms. 

From a diagonal deformation tensor Ji = li+Wi all the 
usual dynamical quantities can then be easily expressed 



-7 



1 



1 



[l + Wi){l + W2){1 + Ws) 



^ a(l + Wi) 



1 



a(l + Wi) 3 



^ i — ^ a — 



Wi 



(1+Wi) 



(12) 
(13) 

(14) 
(15) 



where 9 and a are respectively the velocity divergence 
and the shear tensor. It is important to note that with the 
above expressions a major part of the dynamics is then al- 
ready exactly described. The conservation-evolution equa- 
tions of S (continuity equation), of 9, (Raychaudhuri equa- 
tion) and of (7, are in fact automatically satisfied. We only 
have to find the evolution equations for the Wi compo- 
nents. Such an equation will be derived from an local ap- 
proximation of the non-local Poisson equation. 



Equations (11) give: 

Tr(J-) = -TriJ-'dU) = /3(r)(i - 1) (16) 
In terms of Wi this equation can be written: 

^ {(1 -I- Wj +Wk + WjWk)dlwi 

1=1,3 

= /3(r)(l + ^ + ^)«,.} (17) 

where {i,j,k) is a circular permutation of (1,2,3). This 
equation docs not contain enough information in order 
to evaluate the evolution of each Wi component. How- 
ever its re-writting, respectively, in one-dimensional, two- 
dimensional and three dimensional geometries suggest rea- 
sonable physical approximations, from which we will de- 
duce an approximate evolution of the deformation tensor. 
In a one-dimensional geometry, equation (17) reduces 

to: 



d^wi = I3{t)wi ^ wi{q,t) = a{t)w°{q) 



(18) 



where Wi{q) is the Lagrangian deformation field linearly 
extrapolated to a = 1. This solution is exact for a planar 
collapse. A similar evolution for the three components of 
W tensor corresponds to the Zeldovich approximation. 

In a two-dimensional geometry, equation (17) can be 
written, where (i,j) is a circular permutation of (1,2) 



J2{{^ + Wj)dlwi 



i=l,2 



:/3(r)(l + ^)«;J 



(19) 



In this equation wi and W2 have a symmetric role. Split- 
ting symmetrically this unique equation in two, we get: 



(1 + Wj)dlwi = l3{T){wi + WiWj) + Bi 



(20) 



In these equations Bij components are unknown, however, 
they must obey the following constraints: 

1. Bij = —Bji because the sum of equations (20) must 
reduce to equation (17). 

2. Wj = implies Bij = in order to be consistent with 
the one dimensional case. 

Wi = Wj implies By = 0. In this case we have only one 
unknown and equation (17) gives the solution. 

Making the hypothesis that Bij can be expressed as a 
polynomial of Wi the previous constraints gives: 

— 3. =^ Bij cx {wi — Wj) 

— 2. => Bij oc Wj{wi — Wj) 

— 1. Bij oc WiWj{wi — Wj) 

This means that Bij is at least of third order in Wj. 
Equation (20) with Bij = is then exact up to second 
order included in Wi terms. It can therefore be used to 
compute the evolution of the deformation tensor with a 
very good accuracy. 
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Following the same procedure in three dimensional 
case we get the three following equations (i = 1, 2, 3): 

(1 + Wj +Wk+ WjWk)d\wi 

= + + + -^yfe (21) 

This time the constraints on Byfe are: 

1. Bijk = Bikj , because wj and Wk play symmetrical roles 

towards w, . 

2. Wj = wjc — implies Bijk = in order to be consistent 
with the one dimensional case. 

3. Wi = Wj = Wk implies Bijk = 0, because in this case, 
which corresponds to a spherical geometry, there is 
only one unknown and equation (17) then gives the 
exact solution. 

4. Wk = implies Bijk = Bij in order to be consistent 
with the two dimensional case. 

5. Bijk + Bjki + Bkij — 0, because the sum of the equa- 
tions (21) must reduce to equation (17). 

Making again the hypothesis that Bijk is a polynomial 
in terms of w^, the polynomial of lowest order satisfying 
the constraints 1.) to 5.) is Bijk = B^j + Bik- 

Equations (21) with Bijk = 



(1 + Wj +Wk + WjWk)d^Wi 



/3M(1 



Wj + Wk WjWk 



w,. 



(22) 



are then exact at least up to second order included in Wi 
terms. They form a close system of perturbed equations. 
It is important to notice that our approximation is both 
perturbative and geometrical. It can be considered as per- 
turbative because we discard a term of at least third order 
in Wi in equation (21) and it can also be considered as a 
geometrical approximation because it is exact for planar, 
cylindrical and spherical collapse and only approximate in 
the general case. In the next section we discuss the phys- 
ical meaning of the approximate dynamics of equations 
(22) and we show that our solution reproduces quite faith- 
fully the nonlinear dynamic, in particular the evolution of 
density contrast. 

3. Physical interpretation of the approximation 

3.1. Spherical and ellipsoidal collapses 

The equations (22) which reduces respectively in one- 
dimensional and two-dimensional geometries to equations 
(18) and (20) (with Bij = 0), reproduces the exact planar 
and cylindrical dynamics. In spherical geometry, equations 
(22) also gives the exact solution. As a matter of fact, in 
such a case, there is only one free parameter, the density 
contrast 6 which is directly related to the unique compo- 
nent of the deformation tensor, Wsphericai = ij^Y^"^ ^ 1- 
Such a relation allows one then to identify the collapse of 
a homogeneous sphere with a density contrast 5 and a col- 
lapse in one point with a spherical environment defined by 



Wsphericai whose evolution is given by equations (22). In 

the following, in a similar way, we will identify finite scale 
collapses with local collapses with the same geometry. 

When the three diagonal components of W are differ- 
ent, it is natural to generalize the sphere to an ellipsoid. 
The equations (22) then describe the collapse of a point 
with an elliptic environment which means that the defor- 
mation tensor at such a point can be identified with the 
deformation tensor at the center of an homogeneous el- 
lipsoid. The question is now to find such an ellipsoid. An 
ellipsoid of finite size is defined by three free parameters 
as, for example, the two axis ratios, and the initial density 
contrast (the local dynamical quantities at the center of 
the ellipsoid do not depend on the size of the ellipsoid but 
only on its shape) . These three free parameters are not di- 
rectly related to the initial eigen values of the deformation 
tensor. Consequently in order to link the ellipsoid of finite 
size and the initial deformation tensor we have used the 
tide tensor. This one is very meaningful for the dynamics 
and it is well defined in both cases. The tide tensor at the 
center of the ellipsoid can be computed using elliptic in- 
tegrals, and it is related to W by J" = -(1-|- W)dlW. In 
Figures 1 we compare the evolution of the density contrast 
computed from equations (22) and for a homogeneous el- 
lipsoid having the same initial tide tensor. In principle 
the ellipsoid generates inhomogeneities in the outside ho- 
mogeneous background and these inhomogeneities in re- 
turn perturb the density inside the ellipsoid. The dynam- 
ics of the homogeneous ellipsoid is computed by neglect- 
ing such a feedback efi^ect which means that the density 
inside and outside the ellipsoid remains homogeneous. A 
detailed analysis of the ellipsoidal model can be found in 
V. Icke 1973, White S. & Silk J. 1979. In Figure 1, we see 
that our approximation reproduces very accurately the 
evolution of the density contrast at the beginning of the 
non-linear regime even for very non spherical cases (Figure 
lb). The error is always of a few percent for 5 = 10. In the 
very non-linear phase (for 5 « 100) the dynamical delay 
of our approximation is, in all cases, less than 5%. Figures 
Ic-ld show the final stage of the collapse. The collapse 
time when the density diverges, is estimated with an ac- 
curacy better than 5%. This one remains rather constant, 
whatever the ellipticity. 

This generalization of the spherical collapse allows one 
to describe points with an anisotropic deformation tensor. 
The matter distribution around the fluid elements is ap- 
proximated by an ellispsoid which can fit a generic matter 
distribution much better than a sphere. The picture of an 
ellipsoid is helpful in understanding the dynamics but one 
should not forget that we are evolving Lagrangian fiuid 
elements and not a finite sized ellipsoid. Another impor- 
tant point is that for an elliptic collapse the tide tensor 
remains diagonal, this was our first approximation to ob- 
tain equations (22). Consequently, if we restrict ourself to 
this kind of collapse our only approximation is Bijk = 
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Fig. 1. This figure shows the evolution of the density contrast for two diffcreut ellipsoids. The y-axis is labeled in terms of the 
density contrast and the x-axis in terms of a = a * ^i/a, where a is the scale factor and 5i and a, are the initial density contrast 
and scale factor. The axis ratios of the ellipsoids are 1.24-0.8 for a) and c) and 2.21-0.45 for b), and d) 




Spherical Ellipsoidal Plan-Parallel Parabolic and hyperbolic 

{X='k='^,<0) (^,,X„A,3<0) {l,<0,X='k=0) ( A,,,A,.<OA,>0) (X,<0, X„X3>0) 




Fig. 2. This figure shows the geometrical classification of the different possible collapses, a represents the eigenvalue of the 
shear on the fastest collapsing direction and Oc the collapse time of the fluid element. (The initial density contrast is supposed 
to be positive and equal in all cases). 
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and we have shown that discarding this term was a very 
good approximation. 

3.2. Others Geometries 

The evolution described by equations (22) is not restricted 
to the Lagrangian fluid elements with an elliptical envi- 
ronment for which the three diagonal components of W 
are of the same sign (e.g. the opposite sign of 6) . It can be 
applied to a wider class of initial conditions corresponding 
to an initial W with eigenvalues of different sign. We call 
parabolic (resp. hyperbolic) the points which have an ini- 
tial deformation tensor with two negative (resp. positive) 
and one positive (resp. negative) eigenvalues. For these 
points which correspond to more complex geometries it is 
not possible to find a simple distribution of matter with 
which they could be compared. They might be viewed as 
ellipsoids subject to a strong external tidal fields caused 
by surrounding matter. If the tidal field is strong enough 
to make one direction expand instead of collapsing the 
global matter distribution is not well describe by an ellip- 
soid. The parabolic and hyperbolic geometries allows to 
describe such cases. 

The transition from elliptic to parabolic collapse is the 
planar [i.e. Zeldovich) collapse. As illustrated in Figure 2, 
the shear on the fastest collapsing direction increases and 
the collapse time decreases when one goes from elliptic 
to parabolic and finally to hyperbolic points. Parabolic 
and hyperbolic points collapse even if they are initially 
underdense due to the effect of the shear and of the tide. 
This illustrates that these two quantities, and not only 
the density contrast, play a key role in the dynamics. 
The only points which do not reach infinite density are 
elliptic points which are initially underdense or equiva- 
lently points with three positive eigenvalues for their ini- 
tial deformation tensor. As in the case of the Zeldovich 
approximation, all the points with one initial negative 
eigenvalue for their deformation tensor, collapse. These 
points represent 92% of all points. This ratio can be com- 
puted by integrating the distribution function of the de- 
formation tensor eigenvalues obtained by Doroshkevich 
(Doroshkevicli 1970). Excepting cylindrical and spherical 
points which are highly symmetrical, the direction with 
the most negative initial eigenvalues of W collapses first 
and alone. Consequently the collapsed objects will have a 
pancake shape. 

In the next section we concentrate on elliptic collapses 
and make a detailed comparison of different approxima- 
tions using not only the density contrast but also the shear 
and the tide tensor. 

4. Comparison of different approximations: 

In the previous section we have seen that our approxi- 
mation reproduces with a good accuracy the evolution of 
the density contrast at the center of a homogeneous el- 



lipsoid (for fluid elements with an elliptic geometry). We 

will now push further this analysis in two directions. First, 
we will compare our approximation (hereafter VB as Van- 
ishing Bijk), the one proposed by Bertschinger and Jain 
(Bertschinger & Jain 1994) (hereafter VH as Vanishing 
H-ij) and the first and second order developments (these 
developments are obtained by the method presented in 
the Appendix) to the homogeneous ellipsoid model. Since 
all these approximations keep the deformation tensor di- 
agonal, it is natural to compare them with the dynamics 
of the homogeneous ellipsoid. Secondly, we will show that 
the density contrast alone can be misleading in determin- 
ing the exactitude of a given dynamic. Consequently we 
will also study the evolution of the shear and tide tensor 
which, as we have seen, have a very important dynamical 
role and give valuable information about the geometry of 
the collapse. 

To make the discussion more explicit we describe 
briefly the approximation suggested by Bertschinger and 
Jain. These authors have also proposed a closed local set of 
Lagrangian equations to describe the evolution of a fluid 
element. This set of equations for the irrotationnal case 
and for diagonal shear and tide tensor can be written as: 

5 + a{l + 5)6 = (23) 
e + -e+'^e'^ + aa'oi = -AnGpa^S (24) 

+ -cTj + ^9(Ji + aaf - ^{a'^ak) = -a£i (25) 
a 6 6 

Si - ^f) + -\e£) + Ho-'lfy - ^^a)£) = -A7r^;pH^(7)(26) 

where £ij is the traceless part of the tide tensor J^ij = 
Va;. Vj;^!* {£ij = Tij — ^V^l*) and the dot denote deriva- 
tion with respect to r. Equations (23)-(25) are easily de- 
duced from (7)- (9). Equation (26) governing the evolution 
of the tide tensor is obtained in the Newtonian limit, from 
general relativity by neglecting the magnetic part {Hij ) of 
the Weyl tensor (we will hereafter call the term dis- 
carded in order to obtain equation (26) divided by the sum 
in absolute value of all the terms of this equation). The 
Newtonian expression of the Weyl tensor has been abun- 
dantly discussed recently (Kofman & Pogosyan 1995, 
Bertschinger & Hamilton 1994) but its dynamical influ- 
ence remains unclear. During our analysis we will give a 
partial answer to this question. 

Let us first examine the case of an ellipsoid which col- 
lapses first in one direction, the two other directions being 
equal. The ellipsoid is like a thick disc which collapse along 
its axis of symmetry. The evolution of the density contrast 
is perfectly described both by VB and VH even in the non- 
linear phase {5 w 50) (Fig. 3 a)). The improvement com- 
pared to first and second order developments is significant, 
which illustrates that integrating perturbed equations is 
more accurate than finding a perturbative solution. At the 
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Fig. 3. This figure displays the evolution of different dynamical quantities for the homogeneous ellipsoid model (full line), for 
VB (dashed line), for VH (doted line) and for first and second order perturbative solutions (dashed-doted line). The x-axis is 
always labeled in terms of a (cf. fig 2). Graph a) and b) shows the evolution of the density contrast with diff'erent scales, graph 
c) and d) shows the evolution of the eigenvalues of the shear and tide tensor and graph e) the evolution of Mi as defined in the 
text. The initial axis ratio of the ellipsoid are: 1-0.58. 




Fig. 4. This figure displays tlie same graplis as in figure 3, but for an ellipsoid with an initial axis ratio: 2.4-1 



collapse both VH and VB have a dynamical delay of a few 
percent (Fig 3b)). VH undcr-cstimatcs the collapse time 
and is slightly more accurate than VB which on the con- 
trary over-estimates the collapse time. The analysis of the 
shear components gives similar conclusions. The shear is 
well described at the beginning of the non-linear phase and 
has a small dynamical delay at collapse (Fig 3c)). How- 



ever, the evolution of the tide tensor gives quite a different 
point of view (Fig 3d)). First, wo sec that VH docs not re- 
produce correctly the evolution of £ij even at a very early 
stage. In the linear regime $ oc which implies £ij = 0. 
Thus, the first non-vanishing terms in £ij are at least of 
second order. Consequently the dynamics of Sij is much 
affected by Hij which is generally not vanishing at second 
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order. On the other hand, in VB the term which is ne- 
glected is at least of third order it is therefor normal that 
the correct evolution for is recovered at early stage. In 
the quasi-hnear regime VB still gives quite a good descrip- 
tion for the evolution of the tide tensor while VH does not 
reproduce the correct behavior. All these features can be 
understood by looking at the evolution of Mi (Fig 3 e)) 
Mi grows very rapidly in the linear regime which trans- 
lates immediately in a wrong evolution equation for the 
tide tensor for VH. VB gives the correct evolution of Mj 
in the linear phase. Near the collapse, since one direction 
collapses first, the fluid element evolves towards a planar 
geometry. For this reason Mi is again vanishing (compared 
to other quantities) at the collapse. VB is unable to follow 
this behavior and largely over-estimates the amplitude of 
Mj at the collapse. 

Now we will turn our attention to the case where two 
directions are equal and collapse first. The ellipsoid then 
has a cigar shape and collapses in the two directions per- 
pendicular to its axis of symetry. For this particular ge- 
ometry VB is extremely close to the exact (ellipsoidal) 
dynamics. The density contrast and the shear tensor are 
computed to a high accuracy (Fig 4 a), b) and c)). For 
both these quantities the dynamical delay at the collapse 
is less than 1%. VB is excellent in this case because the ini- 
tial ellipsoid has a cigar shape and its dynamics is close to 
that of a cylinder (VB is exact for a cylindrical collapse). 
Looking now at the tide tensor we see that the two di- 
rections which collapse first are perfectly described while 
the third one, which should be zero for a cylinder, is not 
very well followed. VH does not give accurate results in 
this case. The fluid element tends, at collapse, towards a 
cylindrical geometry for which the Weyl tensor is not van- 
ishing (contrary to what said by Bertschinger and Jain in 
their abstract). Consequently, Mi is playing an increasing 
dynamical role from the linear phase to the collapse which 
makes equation (26) a bad approximation. VB is extremly 
accurate for this particular ellipsoid because its dynamics 
is very close to that of a cynlindar and we have shown 
that VB was exact in cylindrical geometry. 

The two cases that we have studied previously were 
highly symmetrical. They were instructive in order to un- 
derstand the role of the different dynamical quantities in 
each of the approximations. However, they are not rep- 
resentative of generic cosmological fluid elements. Conse- 
quently, we will now study the collapse of an "ordinary" 
ellipsoid with three different axis. We have plotted in fig- 
ure 5 the evolution of only one such ellipsoid but the de- 
scription we give is general and is the result of the analysis 
of many different cases. For these geometries both VH and 
VB describe the evolution of the density field very accu- 
rately (Fig 5 a) and b)). VB is a little more accurate in 
the non-linear phase and significantly better for the es- 
timation of the collapse time. The error on the collapse 
time obtained from VB is always less than 5%, while it 
can be much greater for VH. VB always overestimates the 



collapse time while it is underestimated by VH. The im- 
provement compared to the second order perturbative so- 
lution is in all cases very important. This can be explained 
by the fact that a perturbative solution of equation (21) 
only preserves the perturbative aspect of our approxima- 
tion but loses its geometrical properties which are "coded" 
in the equations. A perturbative solution is inexact for the 
cylindrical and the spherical collapse. The quality of the 
description of the shear and of the tide given by VB is 
quite similar to that of the density contrast (Fig 5c) and 
5d)) . Both these quantities are very well followed at the 
beginning of the non-linear phase and have a small dy- 
namical delay at the collapse (this delay is a little greater 
for the tide). Eventhough our approximation reproduces 
quite faithfully the dynamical evolution of all the usual dy- 
namical quantities, it is less accurate for Mj. In the first 
phase of the collapse the three components Mi start to 
grow, this phase is very well followed by VB. However, in 
a second phase, the Mj components return to zero under 
the correct dynamics (ellipsoidal model) while they con- 
tinue to grow in VB (Fig 5 e)). As was mentioned before, 
when one direction collapses first the ellipsoid approaches 
a pancake shape when reaching the collapse and this ex- 
plains why Mj evolve towards zero. However, even if they 
tend to zero at the collapse the Mj terms have a very im- 
portant dynamical role. Eventhough Mj terms are small 
(but not vanishing) in the linear regime, they increase very 
rapidly to become one of the dominant terms in the equa- 
tions (26). Consequently, VH rapidly becomes a very bad 
approximation. This explains the rather strange behavior 
of the tide tensor predicted by VH (Fig 5 d)). The error 
initially present in equation (26) propagates to equation 
(25) and then to (24) and finally to (23). One component 
of the shear is much affected by this error (Fig 5c)), VH 
predicts two negative components for the shear while there 
is only one. The density contrast is less sensitive to it. 

From the analysis of the different approximations we 
can conclude that the Mj terms slow down the collapse 
and couple the directions. VB which overestimates the 
amplitude of Mj always collapses too slowly and the three 
components of its tide tensor are much closer than in the 
exact ellipsoidal dynamics (Fig 5 d)). On the contrary, 
VH which completely discards Mj, collapses too quickly. 
Figures 5 c) and d) illustrate the importance of Mj for the 
dynamics and for the geometry of the collapsed objects. 

5. Conclusion 

In this paper we have presented and tested a new lo- 
cal Lagrangian approximation to the dynamics of a self- 
gravitating cosmological fluid. One originality of our ap- 
proach was to describe the gravitational dynamics through 
only one quantity, the deformation tensor. We have shown 
that with this tensor it is possible to write the dynamical 
equations in a synthetic form and to highlight the role of 
different important dynamical quantities, especially the 




Fig. 5. This figure displays tlie same graplis as in figure 3, but for an ellipsoid with an initial axis ratio: 1.24-0.8 



shear and the tide. The evolution of the deformation is 

evaluated by deriving a close set of local crjuations from 
the exact dynamical system after two approximations. The 
first one supposes that the deformation tensor remains 
diagonal during its evolution, the second one consists in 
discarding terms of third order in terms of the deforma- 
tion tensor component. The approximate dynamics thus 



obtained, which reproduce exactly the planar, cylindrical 
and spherical collapses is then exactly integrated. Such a 
procedure is completely different from the method which 
consists in deriving a perturbative solution from the exact 
dynamics. 

In order to test our approximation we have compared 
it with the collapse of a homogeneous ellipsoid and with 
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the approximation proposed by Bertschinger & Jain. Prom 

this analysis wc conclude that our approximate dynam- 
ics is a very accurate tool for following all the dynami- 
cal quantities during the collapse. We have also, for the 
first time, shed some light on the important dynamical 
role played by the Newtonian counterpart of the magnetic 
part of the Weyl tensor, discarded by Bertschinger & Jain. 
This tensor slows down the collapse and couples the di- 
rections, discarding it reduces the collapse time signifi- 
cantly and greatly modifies the geometry of the collapse. 
The main weakness of our approximation is that, unlike 
the Zeldovich approximation, it is unable to reconstruct 
the trajectories. However many cosmologically interesting 
quantities, such as for example the mass function, can be 
computed without reconstructing the final Eulerian field. 
Furthermore one of the main interests of this work might 
is that it may give a better understanding of the very 
complex phenomena which occur in the non-linear phase 
of the gravitational dynamics. 

We would like to thanks M. Maccormick for rereading 
the manuscript. 

6. Appendix: Link with perturbative theory. 

In section 4 we have seen that integrating our system is 
more accurate than finding a perturbative solution, at first 
and second order, of equations (22). Would this still be 
true for higher order solutions and does the perturbative 
solution of n*^ order converge toward the exactly inte- 
grated solution of equation (22)? In order to answer these 
questions we have developed an analytical recursive for- 
mula which relates the perturbative solution of equation 
(22) at order n to the solution at order (n — 1). The con- 
vergence of this development and the ability of such a 
perturbative solution to reproduce all the features of the 
exactly integrated solution will be the object of a forth- 
coming paper. 

In this section we give a recursive formula to compute 
the solution of our system at any order. Of course this 
development will differ from the development of the ex- 
act dynamics beyond second order, but by going to higher 
orders we can analytically recover the accuracy of the nu- 
merically integrated solution. 

Equations (21) can be written, after some simple alge- 
bra, as: 

(1 -I- Wj +Wk+ WjWk){d\wi - f3{T)wi) 

= -f3ir)i^^^^ + lwjw,)w, (27) 

It is possible to find a perturbative solution at any order 
of the preceding equation. Let us first introduce a few 
useful notations: 

oo 

E(n) 
wl 

n=l 



n-1 

(n) (;) (n-l) 

1=1 

n — 2 n — m — 1 
a(n) \ ^ \ ^ (m) (I) (n-l-m) 

m=l 1=1 

Pi(") j2 (n) Qi \ in) 

With this notation it is possible to rewrite equation (27) 
as: 

oo n— oo 

n=l n=l 

CO {■!'■) I (r?) 

= -/3(r)E( 2"" +3^'"'^ ^^^^ 

n=l 

Keeping only terms of n*'' order in this equation gives: 

pin) ^ _^^^^(-^t±^Apin)^^Y^^^^^ 
^ 1=1 

(n) (k) 

Pj can be computed using only the with k < 

n-l. The differential equations d^Wj'"'' - P{t)w^"^ = p/"'' 
can be recursively integrated to any order. The fact that 
p/""* is of order (n) (except Pj^"^ which is zero) is also 
shown recursively. If we keep only the linear growing mode 
at first order, these equations can be integrated very easily 
to a high order. Since equation (27) is exact only at second 
order the development will not be exact beyond that order. 

However, as we have seen, the exact solution of equa- 
tion (27) is much more accurate than the second or- 
der solution and by developing to a higher order we can 
have an analytical solution very close to the numerically- 
integrated one. Another interesting feature of this devel- 
opment is that since we can go to any order very easily it 
is possible to check if a given quantity can be accurately 
approximated by perturbative theory. It is possible that 
some quantities are not well described by perturbative the- 
ory and are only computable with the exact solution. For 
example d^Wj becomes infinite at collapse for the exact 
solution but any perturbative calculation will give a finite 
value. 
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